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Abstract 

The exciton dynamics in a molecular heterodimer is studied as a function of differences in excitation and reorgani- 
zation energies, asymmetry in transition dipole moments and excited state lifetimes. The heterodimer is composed 
of two molecules modeled as two-level systems coupled by the resonance interaction. The system-bath coupling is 
taken into account as a modulating factor of the energy gap of the molecular excitation, while the relaxation to 
the ground state is treated phenomenologically. Comparison of the description of the excitation dynamics modeled 
using either the Redfield equations (secular and full forms) or the Hierarchical quantum master equation (HQME) 
is demonstrated and discussed. Possible role of the dimer as an excitation quenching center in photosynthesis self- 
regulation is discussed. It is concluded that the system-bath interaction rather than the excitonic effect determines 
the excitation quenching ability of such a dimer. 
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1. Introduction 

Photosynthetic pigment-protein complexes arranged as light-harvesting antenna absorb the Sun light and trans- 
fer the captured energy towards the photosynthetic reaction center, where electron transfer across the membrane 
is initiated. In addition light harvesting antennae also carry out self-regulation function, which is a physiologically 
highly significant strategy evolved by plants [1 . Thanks to the excitation energy density control in photosystem II, 
termed non-photochemical quenching (NPQ), plant photosynthesis can function efficiently under very different light 
conditions, from low up to very high intensities corresponding to daily changes of illumination at the same location. 
The NPQ phenomenon is usually attributed to some activated quenching species which allow the excitation to 
undergo a rapid non-radiative decay. The exact location of the quencher within the antenna and its precise nature 
are a matter of on-going debate, with both chlorophylls (Chi) and carotenoids (Car) being put forward as essential 
components of the quenching mechanism. The promising candidates for the quenching mechanisms (in no order 
of preference) are: (i) quenching by a Chl-Chl dimer showing the charge transfer (CT) state character |5] or via 
the formation of Chl-Chl excimeric states [3^ ; (ii) by a CT state between Chi and Car (xanthophyll) resulting in 
generation of the cation radical state of the xanthophyll [4l|5]; (iii) via the excitonic coupling of Chi to a short-lived 
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xanthophyll excited state |3] , and (iv) by the direct energy transfer from the Chi pool to a particular Car (lutein) 
[7] . Both latter suggestions are based on the fact that the lifetime of the lowest excited state of the lutein molecule 
is short (<10 ps). Evidently, this could explain the possible quenching efficiency in the light-harvesting complex 
of photosystem II (LHCII) if a specific arrangement between the Chi and Car molecules is established under the 
NPQ conditions. Thus, according to the latter suggestions, the simplest candidate which might be responsible for 
the NPQ, is an asymmetric molecular dimer (a heterodimer). The possibility to attribute the excitation quenching 
ability to the Chl-Car heterodimer was also suggested and experimentally studied for model dyads composed of 
Car and tetrappyroles |H1 E] . A dimer is also the simplest molecular aggregate where the excitonic features are well 
expressed. 

Spectral characteristics of molecular dimers are often remarkably different from those of the individual molecules. 
Their absorption bands can be significantly shifted in comparison with those of their constituent molecules, mainly 
due to intermolecular interactions resulting in a nonlocal character of their excited eigenstates. These delocalized 
eigenstates of electronic excitations are usually termed molecular excitons [TUl fTTI fT^ [T^ . The interaction of 
electronic excitations with intra- and inter-molecular vibrations causes a disruption of the phase relationship between 
excited states of the molecules constituting the exciton wave functions |13l ITl]. Such type of interaction makes a 
distinct influence on the exciton dynamics, and plays the dominant role in determining the exciton relaxation 
pathways. 

Usually, a homodimer is used for considering various aspects of the exciton dynamics and relaxation |15l I16L 
[171 [18]. For a heterodimer, the distinctness of constituent monomers is often limited to excitation energies |19} [20]. 
Some aspects which could be attributed to the heterodimer were also disclosed by analyzing the exciton-CT state 
mixing problems |21l [55]. Recent experiments based on coherent photon echo measurements demonstrated the 
possibility to follow the coherent phase dynamics and incoherent population relaxation of excitons in photosynthetic 
pigment-protein complexes [TOl 1551 E^ . However, possible role of the coherent dynamics in determining such type 
of the regulation function - attributing it to the ability of the excitation quenching in heterodimers - has not 
been considered yet. Here, such type of the analysis is presented. We consider the effects originated from the 
differences in the excitation energies, reorganization energies and excitation lifetimes of the constituent molecules 
of the heterodimer. 

The paper is organized as follows. In section [2.1| we introduce a theoretical model of an excitonic dimer and 



discuss the general form of equations describing the excitation dynamics in the system. In sections |2.2| and 2.3 we 
introduce two specific theories of relaxation that are nowadays standard in the literature and in section ]2T4] we discuss 
our model for the radiationless excitation energy dissipation. Numerical results for the excited state dynamics under 
various parameters of the dimer system with and without relaxation to the ground state are presented in sections 
3.1| and |3.2[ respectively. The consequences of these results and the possible role of the heterodimer in determining 
the NPQ are discussed in section |4| 



2. Modeling 

2.1. Excitonic Dimer 

The exciton spectra of molecular aggregates are usually considered in the so-called Heitler-London approximation 
[lOl [TTJ [T^. According to this theoretical concept, monomers comprising the dimer are characterized by two 
electronic levels, which reflect the resonant transition in the spectral region under consideration. The states of the 
dimer can be constructed on a Hilbert space, which is a direct product of the Hilbert spaces of the monomers. The 
ground state of the dimer is constructed from monomers i = a and b as |g) = \ga)\9b)- Similarly, two singly-excited 
states are deflned as \i) = \ei)\gj), j ^ i where the ith molecule is in the excited state and the jth molecule in 
the ground state. The two singly-excited states are coupled to each other by the resonance interaction Jab, and 
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consequently they are not eigenstates of the system. A double-excited state |/) = |ei)|ej) where both molecules are 
excited is not relevant in our analysis. The given states form the site basis. 

Interaction of the excitonic system with its environment forming a thermodynamic bath, has to be taken into 
account when analyzing the exciton dynamics and relaxation. Thus, the full Hamiltonian H should contain a purely 
system part Hs, describing only the electronic excitations, the bath part Hb and the system-bath interaction part 

H = Hs + HsB + Hb. (1) 
The individual terms in the site basis explicitly read: 

Hs = + ^^)^' + Jab{\a){b\ + \b){a\), (2) 

i 

HB=Tip) + Vg{q), (3) 
iJss -^AV^,(g)if,. (4) 

i 

Here, i = a, b; e° denotes the energies corresponding to the electronic excitation in the zth monomer, and Ki = 
is the projector onto the zth site. The system Hamiltonian Hs has the form of the Frenkel exciton Hamiltonian. 
Operators T{p) and Vi{q) denote the kinetic energy of the nuclei and the nuclear potential energy surface of the 
ith site, respectively; p and q denote the generalized momenta and coordinates of the bath. We defined the 
reorganization energy as 

A, = mq) - Vg{q)),, (5) 
where the angular brackets (. . .)q represent the averaging over the equilibrium bath: 



{Re,..}- (6) 

Here, Trq {. . .} denotes the trace operation, and R^q is the canonical equilibrium density operator of the bath: 

R,q^e-P"-/Trq{e-^"-}, (7) 

where /3 = \/k^T, T is the temperature and is the Boltzmann constant. The so-called energy gap operator 
/S.Vi{q) defined as 

AVM)^V,{q)-Vq{q)~\ (8) 

describes the fluctuations of the energy gap between the potential energy surfaces of the state \i) and the ground 
state \g). We next denote = -I- A^. 

In order to find the excitonic states, which diagonalize the system Hamiltonian, we use a unitary transformation: 

i 

which yields energies of the exciton states: = (U^^ HsU) In the following we use the Roman letters for the 
site basis and the Greek letters for the eigenstate basis. In the case of a dimer the transformation matrix has a 
convenient analytical form [2S] : 
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\^ sm^ cost' ) ^ ' 

where d is the so-called mixing angle, which is defined as 

6 = - arctan — (11) 

and has values ranging from to 7r/2. 

The most detailed description of our system of interest and its environment is given by a density operator W 
acting on the joint Hilbert space of system and bath states. However, since we are interested only in the electronic 
degrees of freedom (DOF), we switch to the reduced density operator (RDO) /5, which is obtained by tracing the 

full density operator over the irrelevant (environment) DOF as p = Tvq the basis of electronic states 

the RDO is represented by its matrix elements pab = {o\pW)- Usually it is assumed that the initial full density 
operator Wq can be factorized into the RDO, po, and the complementary equilibrium density operator of the bath: 
Wq = pq (g) Rf.q. This factorization assumption is well motivated in the spectroscopic applications by the existence 
of a single electronic ground state and a large energy gap between the ground and the excited states. Ultrafast 
excitation of the system by laser light corresponds to the Franck-Condon transition (i.e. the nuclear DOF are not 
involved in the process), and the state remains factorized immediately after excitation. 

Using these definitions the equations of motion for the RDO are given in a general form by [TT, "25] : 

^^p=-iCsP + V[p]+K.p, (12) 

where on the right hand side we have three superoperators acting on the RDO. The first of those is the Liouville 
superoperator associated with the corresponding system Hamiltonian (Eq. (|2])) by the relation Cs» = \ [i?5,»]. 
The square brackets denote a commutator as usual, and h is the reduced Plank's constant (from now on we set 
it to unity for convenience). This part of the equation governs the coherent evolution of an isolated system of 
electronic states. The second term represents the general propagation scheme for the electronic system due to the 
interaction with the bath DOF as detailed below. The last term is introduced to treat the decay of excitons due to 
the nonradiative transitions from the excited states to the ground state. 

2.2. Redfield Relaxation Scheme 

One of the most popular methods of treating an excitonic system coupled to a bath are the Redfield equations. 
They follow from the Liouville-von Neumann equation under the assumption that the system-bath interaction is 
sufficiently weak to perform a perturbative expansion. In the second order approximation in the interaction Hamil- 
tonian HsB one gets the so-called Quantum Master Equation [14 . Then, performing the Markov approximation in 
the excitonic basis we get the dissipative part of Eq. reading: 



(-D[p])pi. = -^7^^^.,^'^'PpV'■ (13) 

The tetradic relaxation matrix TZfj^u.u'v'it) is the Redfield tensor and it is most conveniently given as follows: 



72-/ji/,/jV (t) —Syu' ^ ^ r^g.j^' (t) -|- Sfip^i ^ ^ ^ue.,ev' (^) 



(14) 
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Here, 6^,^' is the Kronecker's delta, and F's are certain damping matrices defined as follows: 

T^,,^,,,{t) = y2{^l\K,\l^){^^'\KJ\l,') f drC,,(r)e''^'^V-. (15) 

The matrix elements {filKilv) represent the basis transformation from molecular states (sites) to the delocalized 
eigenstates, and Cij{t) is the energy gap correlation function in the site basis 

a,,{t) ^ {^vdt)^vM),- (16) 

The operator /Syi{t) — U^{t)AViUB{t), where UB{t) is the bath evolution operator, represents the evolution of 
the energy gap driven by the bath Hamiltonian Hb- In the following we make an assumption that the energy gap 
fluctuations at different sites are uncorrelated, i.e. Cij{t) = SijCi{t). 



Equations (12 1 can be solved once we have an explicit form of the energy gap correlation functions (Eq. (16l) 



C(t) or spectral densities C"(w). The two quantities are related by 

C,it)= / dw{(l + n(t^))e-'"*+n(w)e"^*}C,"(w), (17) 
Jo 

where n{uj) denotes the Bose-Einstein distribution function: 

«(^)^ /\ 1 ■ (18) 

exp(cL>p) — 1 

A well explored model of the bath corresponding to a continuous spectrum of the phonons is the so-called Debye 
spectral density ^ ,264 22^ ■ 

Cr(a.) = 2A,^^, (19) 

where the parameter 7 characterizes the timescale of the dissipation of the reorganization energy, and the relation 
Jp duiC" (ui) / WK = Xi holds. Hence we can characterize each monomer in our system by the reorganization energy 
Xi deflned in Eq. (|5|, which is not just the energy shift in Eq. but also the measure of the system-bath 

coupling. The parameter 7, which is purely the bath property, is identical for both monomers under the assumption 
of identical baths. 

At this stage the Redfield equations can already be used in calculations, however, under the conditions justifying 
the Markov approximation the related electronic dynamics is much slower then the decay of the bath correlation and 



the scheme can be simplified by shifting the integration limit in Eq. ( 15 1 to i — cx). This yields a time-independent 
Redfield tensor. The so-called secular approximation |14| . which decouples the evolution of populations from that 
of the coherences, also applies under these conditions. Formally, this approximation is realized by setting to zero 
elements of 'R-^v,^l>u' with the indices other than those satisfying the conditions ^ = v and fi' = v' or fi = fi' and 
V = v' . In this paper both versions of the Redfield theory - the full and the time-independent secular - are used for 
comparison. 



2.3. HEOM Relaxation Scheme 

Recently the hierarchical equations of motion approach has been introduced for the description of the exciton 
dynamics (see, for instance, [28'). Since it employs the cumulant expansion technique, it is exact if the bath is 
Gaussian. This is the case for a bath of harmonic oscillators. Unlike with the other methods, such as the Redfield 
equations, the precise form of HEOM depends on the specific form of the bath correlation function. Here we use 
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an approximate HEOM theory termed the hierarchical quantum master equation (HQME) |^ . It is based on the 



following form of the correlation function, which is an approximate expansion of Eq. ( 17 1 with the Debye spectral 



density defined by Eq. (19 1 



at) = ("^ - ^) e-* - iAae-* + f6 (t) ; (20) 

where S {t) is the Dirac delta function. 

In this method we replace the RDO p in Eq. ( |12[ ) with a hierarchy of operators p ^ Pn indexed by n = 
(ni,n2, . . . tTIm), where n^'s are non-negative integers and M is the number of states in the system. In this set, 
the operator denoted by po, where = (0, 0, ... , 0), is the RDO. Other operators are auxiliary and they contain 
information about the system-bath correlations. The dissipation term reads within the HQME as: 

M M 



1=1 i=l 
M M 



i=l i=l 

where we have defined the following superoperators: 
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and Bi* = [Ki, •]. The curly brackets denote the anticommutator, and = (ni, 712, . . . , ± 1, . . . , um)- The 



operators pn with certain n such that N = J2i '^ii define a "tier N" of operators. We can see that within Eq. (21 1 
each tier N becomes coupled to neighboring tiers ± 1. Formally, the hierarchy in HQME continues to infinity. 
In practice, however, one has to truncate the hierarchy by setting all auxiliary density operators from the tiers 
N > Ntrunc to zero. The cutoff Ntmnc is defined by the relative values of the model parameters. In calculations it 
can simply be chosen to ensure the convergence of the solution for po ■ 

2.4^. Exciton Decay 



In contrast to the second term in Eq. ( 12 1, which is described by relaxation theories of previous subsections, we 
include decay of excitons in a phenomenological way by taking the experimentally estimated rates for population 
decay processes into account. Therefore, the /C superoperator elements are phenomenologically defined in the site 
basis. We assume that the only non-zero elements of this superoperator are those that describe the decay of 
populations and the corresponding decay of coherences, i.e., we assume the secular approximation. The elements 
of the superoperator thus read: 

K^ijM = — -SikSji, (22) 

where Ki denotes the relaxation rate of the population of the ith site. This way we have only two types of elements: 
the population relaxation (ICuAi) and the coherence decay due to corresponding population relaxation [ICij^ij). 
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absorption 



Figure 1: Schematic representation of the constituent monomers of the heterodimer. The excited states of the monomers are indicated as 
Ea and ei,, correspondingly, the resonance interaction between the monomers is indicated as J^b- The optical transition in the a-monomer 
is allowed, while it is forbidden in the b-monomer as indicated by the vertical black arrows. The inset shows a hypothetical absorption 
spectrum of such system. The position of e;, might be above or below ea and therefore is shown by dotted line; reorganization energies 
are shown to correspond to the homogeneous linewidths of the absorption spectrum. 

In order to use the /C tensor with the Redfield equations, it needs to be transformed into the excitonic basis. 
The transformation matrix for a superoperator is given in the Appendix. In the case of a dimer we obtain simple 
expressions with the mixing angle once again. For the population relaxation the tensor elements read: 

1^130,013 =ICaa,aa sin^ 9 + K-bbM COS^ 9. (23) 

3. Results 

We set up our heterodimer to have some typical properties of photosynthetic Chl-Car aggregates. As discussed 
above, such a dimer might be relevant to the quenching in photosynthesis. Therefore, one of the constituent 
monomers is characterized by an optically dark and extremely short-lived excited state. The other monomer has 
a strong transition dipole moment and a long-lived excited state. Their homogeneous broadenings are different as 
well. This situation is schematically depicted in Fig. [T] 

The constant parameters are the following: T — 300 K (which corresponds to ksT ~ 208 cm~^), Jab — 100 cto~^, 
— 100 /s, dipole moments dag = 1, dbg = in accord with Fig. nl The tunable parameters are the reorganization 
energies Aq, Af,, given as four combinations of the values 30 cm and 150 cm~^, and the gap between the site 
excitation energies with values zblOOcm"-'^ ("-I-" corresponds to the bright state being above the dark one; " — " 
corresponds to the opposite situation). 

We will consider optical excitation as a trigger of the relaxation dynamics. The initial condition for the evolution 
of the RDO corresponding to a resonance excitation of the optically allowed state is then given by p^i_i{Q) = M/igP, 
where d^g is the matrix element of the transition dipole moment operator in the excitonic basis. 

The main parameter characterizing a heterodimer is the difference of site energies, but as can be seen in Fig. 
[2j upon the presence of the bath, in the case of A^ ^ A^, the definition of the energy gap becomes somewhat 
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Figure 2: Definition of the energy gap in the case of different reorganization energies. The dashed vertical arrows indicate the Franck- 
Condon transitions. Clearly Aa ^ \b corresponds to Ae" ^ Ae. 

ambiguous: we can define it either as Ae ~ Ca — ^b or as Ae" = — (even though the relation Ae = Ae" + Aa — 
holds). Application of these different definitions of a dimeric energy gap has a distinct infiuence on the modeled 
dynamics. 

3.1. Population Kinetics 

We first study the dynamics of the excited state in the absence of relaxation to the ground state, i.e. when 
/C = 0. In Fig. [3] we show the time evolution of the higher excitonic state population modeled by the Redfield 
equations (a and b), the HQME (c and d) and the secular Redfield equations (e and f) using four combinations of 
the reorganization energies. The left column corresponds to Ae = +100cto~^ as a fixed energy gap, while the right 
column corresponds to Ae" = +100 crTi^^ being fixed. The equal reorganization energies (black and blue curves) 
serve as a good starting point for comparison of the methods. As we can see, the initial stages of the Redfield 
and HQME solutions (^ 100 fs) look very similar while the rates of the thermalization and the frequencies of the 
coherent oscillations are different (the coherent oscillations are more dramatic but shorter-lived in the Redfield 
solution). Moreover, different methods give us different equilibrium population values, and within the HQME 
solutions the latter do not coincide for two pairs of identical A's. 

The case of different reorganization energies (red and green curves) is quite nontrivial [22|: the results for both 
energy gap definitions obtained by all three methods are considerably different in the long-time limit. In the right- 
hand column the initial values are scattered simply because fixing Ae" with different A's gives us different Ae used 
in the definition of the excitonic basis. Moreover, the red curve in Fig. [3^ reveals the well-known problem of the 
Redfield scheme, namely, that it does not guarantee positivity. The secular Redfield solutions do not suffer from 
this problem, and in this case all four A combinations give relaxation to the same equilibrium position. 

So far we have treated the excitonic states defined by Eq. ^ as the eigenstates of the system. However, HQME 
solutions in the excitonic basis, as given above, exhibit non-vanishing steady state coherences. The presence of the 
steady state coherences can be interpreted as a refiection of renormalization of the system eigenbasis [25^ taking 
place in the course of time. By virtue of the exact treatment of the system-bath interaction the HQME solutions 
are basis-independent, which allows us to identify the so-called preferred basis (or "global basis" [30J in which the 
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Figure 3: Evolution of the higher excitonic state population. The insets show the combination of reorganization energies given in the 
form Aa/Af, (cm~^ /cm~^). The left column corresponds to Ae = +100 cm~^ , the right one - to Ae" = +100 cm~^ . a) and b), c) and 
d), e) and f) correspond to the full Redfield scheme, the HQME and the secular Redfield equations, respectively. 
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t igure 4: The HQME solution of the Ae" = 100 cm configuration in the preferred basis. The inset shows the combination of 
reorganization energies given in the form Aa/Aj, (cm~^ /cm~^). 



density matrix is diagonal in the long-time limit (see, e.g. |31]V This is done by diagonalizing the stationary part 
of the solution in the excitonic basis. Representing the RDO in the preferred basis for the fixed Ae configuration 
has the effect of slightly shifting the equilibrium population values. Their order, however, remains the same as in 
the excitonic basis. Meanwhile for the fixed Ae*^ configuration (Fig. [4|, we have not just quantitative but also 
qualitative changes in comparison with Fig. [3}i. Evidently, the solutions for Ai ^ A2 configurations now converge 
into the same steady state value which is in-between the ones for equal reorganization energies. 



3.2. Relaxation to the Ground State 

We next include the relaxation to the ground state. The relaxation tensor of the equation of motion for the 
RDO is given by Eq. ( 22 1 . We are interested in the influence of the resonance coupling on the lifetimes of the dimer 
eigenstates — n^^. For our heterodimer we assume that in the absence of coupling the excitation of one state 
relaxes considerably faster than the other (rj, <C or JCaa.aa ^ J^bb.bb)- Therefore we can simplify Eqs. (23) to 



'Caa.QQ =K.bbM sin^ 9; (24) 

^P0,I3I3 —JCbb,bbCOS^ 9. 

Denoting the lifetimes of the short-lived state by ts = Tp = and the long-lived one by tl = ^ ^/l^aa,aa 

we present their dependence on the (site) energy gap normalized to the resonance coupling in Fig. |5] The lifetimes 
are normalized to the initial lifetime of the (original, i.e. uncoupled) short-lived state T(,. 

Now let us consider the interplay between processes of thermalization as described in the previous subsection 
and relaxation to the ground state. The numerical results were obtained from the HQME with the relaxation tensor 
and are presented in the preferred basis. The results for a system with Ae — 100 cm~^ and — 2 ns, ti, — 5ps 
are shown in Fig. [6] (the corresponding "excitonic" life-times yield « 21 ps and « 7ps). The populations are 
given in the logarithmic scale in order to reveal the two-exponential nature of the process. Only the higher state 
population evolutions are shown, because the lower state evolves with identical rates. The result is similar to the 
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Figure 5: Excitonic mixing of lifetimes. Indices "S" and "L" denote the short-lived and the long-lived states accordingly. The energy 
gap Ae is normalized to the resonance coupling, and the lifetimes are normalized to the lifetime of the short-lived state in the absence 
of resonance coupling (r;,). 



one shown in Fig. |4] except that instead of the steady state values we have now a decay due to relaxation to the 
ground state, and the long-time population values for the Xa ^ Af, configurations (red and green curves) are no 
longer identical. Obviously, the long-time decay rates are identical in all four cases. 

To demonstrate the influence of the energetic position of the short-lived state, the evolutions for the case 
of Ae — ±100 cm^^ and = Af, = 150cto~^ are considered. For demonstration, we provide an analytical 
approximation of the evolution while neglecting the coherent oscillations. This can be achieved by using a simple 
system of rate equations: 

{X — ^{Kx^kx^y)x + ky^xy\ (25) 

Here, x and y are the populations (of the accordingly lower and higher states), Kx and Hy are the population 
relaxation rates from the tensor /C, and kx^y, ky^x are the effective thermalization rates obtained from the data 
fit of the solution of the HQME without relaxation to the ground state. Following this approach, the population of 
the higher state reads: 

y{t) = Cyie-^'' + Cy2e-^'\ (26) 



where ^i, ^2 are the eigenvalues of Eq. (25l, and C's are the corresponding amplitudes. When kx^y, ky^x 3> Kx, Hy 



the second-order small terms can be neglected thus giving: 




(27) 



This allows us to estimate the dynamics in the following way. First, we determine quantities Cyi, Cy2 and ^2 fi'om 
the HQME solution without K, (hence ^1 — 0) by numerical fitting (disregarding the oscillations). Then we derive 
the rates kx^y and ky^x from the known value of ^2 and the detailed balance condition. It is noteworthy that 
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one must use an effective energy gap instead of the purely excitonic one if the reorganization energies are large. 
Finally, we calculate Kx and Ky from the relaxation superoperator /C in the preferred basis, and then construct the 
^1 eigenvalue. 

The results obtained are depicted in Fig. [7] using the logarithmic scale. The populations of the long-lived state 
(in black) are shown together with their analytical approximations (in green). In the case of Ae = — 100 cm~^ 
(dashed curves) the lower state is the long-lived one (hence, originating from the bright electronic state), therefore 
the thermalization part gives a minor effect. As can be seen, the actual rate of the relaxation to the ground state, 
T, depends on which of the two states - the short-lived or the long-lived - is the lower one in the dimer. 

The analysis given above provides an explanation of the long-time decay rate (^i) sensitivity to the position of 
the short-living state. In particular, swapping the states results in swapping the values of and Ky in Eqs. (25). 
Therefore, the second term in the first equation from Eq. (27 1 changes the sign, and that is the reason why the 
overall decay rate is sensitive to which state - the fast decaying or the slowly decaying - is the lower one. 



4. Discussion 

As follows from our calculations, the asymmetry in reorganization energies allows for two non-equivalent def- 
initions of the energy gap in the one-exciton manifold. These two cases clearly have different physical meanings. 
Since the energy e refers to the Franck-Condon transition region of the potential energy surface (see Fig. [2]), Ae 
corresponds to the distance between the peak positions in the absorption spectra. Therefore we can loosely call it 
the "optical energy difference". This has a perfectly clear meaning in the absorption spectroscopy, however, the 
thermal equilibrium does not establish itself with respect to this energy gap. Should the resonance coupling be 
simply perturbative, we could expect the thermalization with respect to Ae*' as the first approximation. Since this 
gap is usually used in the Forster resonance energy transfer (FRET) theory |13l [T^ we shall call it the "Forster 
energy difference". 

For the fixed optical energy difference, in the case of Ao = Ab an equal amount of energy associated with each 
state is dissipated during the thermalization, and therefore the same equilibrium values are reached. In the case of 
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Figure 7: Solutions of the HQME with relaxation and approximations. The populations are given in the logarithmic scale. Indices "S" 
and "L" denote the short-lived and the long-lived states accordingly. 



Aa 7^ Af), the vibrational relaxation shifts the initial energies by a different amount (Ae" Ae). Hence, according 
to the Redfield relaxation, in Fig. [3^ the energy gap between the states is effectively decreased (green curve) or 
increased (red curve) (here, the increase of the energy gap is seemingly the reason for the failure of the Redfield 
scheme). The case of fixed Forster energy gap (Fig. ^jp) demonstrates the excitonic mixing of reorganization 
energies even better in the sense that for uncoupled monomers the vibrational relaxation would yield identical Ae*^ 
regardless of the combination of A's. The actual situation is schematically shown in Fig. [Sj An interesting situation 
arises in the case of Aa < Af, (Fig. ^jp, red curve), because the two states are swapped in comparison with the other 
combinations (both initially and in the long-time limit). The way all the equilibrium values are ordered (which 
corresponds to the width of the energy gap that determines the equilibration) tells us the peculiarity of the Redfield 
scheme, namely, that the reorganization energies are mixed according to the "initial", i.e. optical energy gap (cf. 
Eq. (111). Hence, the initial excitonic configuration is maintained throughout the whole evolution. 

This is even more strongly pronounced in the secular Redfield scheme (Fig. |3^ and f) where the equilibration 
is fully determined by the optical energy gap Ae alone. This shows us how severely the neglect of the interplay 
between the populations and coherences can distort the picture of thermalization. For instance, it would suggest 
that dimers with identical absorption maxima positions but different line-widths would have identical equilibrium 
population distributions. The full Redfield equations, on the contrary, demonstrate that the optical energy gap, as 
the most intuitive input parameter, along with the resonance coupling are not enough for the description of thermal 
equilibrium in the case of different reorganization energies. Moreover, even the Forster energy gap can be a poor 
benchmark for assessing the equilibrium populations (cf. Fig. [s] bottom panel). 

The HQME solutions (Fig. [sJ; and d) bear some resemblance to the results that follow from the Redfield 
equations. However, an unexpected feature is the difference between equilibrium population values for the case 
when Aa = Ab, which again tells us about different energy gaps in each case. This points to the action of the bath 
upon the system. And indeed we find that the conventional excitonic basis is just an approximate eigenbasis due 
to the non-vanishing coherences. Establishing the preferred basis as given previously led us to conclude that a 
good measure of the action of the bath is an effective resonance coupling J^f / , which has been used to replace the 
original J value with to reproduce the effective energy gap defined by the thermal equilibrium in the eigenbasis. It 
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Figure 8: Schematic representation of the Forster energy gap in the presence of the resonance coupling. Monomeric units with identical 
Ae'^ are shown on the left. The dotted arrows on the right represent the vibrational relaxation (and tilde is used to denote the excitonic 
basis). The top row corresponds to Xa = Xb, the middle and the bottom ones - to Aa > Aj, and Aa < A;, accordingly. 
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appears (cf. Fig. |4| that the connection between the original resonance coupHng and the effective one resembles 
the so-called small polaron transformation [32]. The details of such dynamic suppression of the resonance coupling 
captured by the HQME is beyond the scope of this paper and are presented elsewhere [30] . 

Discussing the relaxation to the ground state requires a couple of remarks about the excitonic mixing of lifetimes. 
From Eqs. ( 23 1 or ( 24 ) we notice that the mixing of lifetimes does not depend on which of the two states is the lower 
one. In the special case of one state being extremely long-lived, we can obtain the following picture: for a dimer 
degenerate in energies both excitonic states would have lifetimes twice the lifetime of the short-lived state in the site 
representation, whereas for different site energies in the limit of negligible coupling J we retrieve the original values 
of the lifetimes. Now, the intermediate regime (energy gap not too big and/or rather strong resonance coupling) 
tells us that the shorter lifetime varies rather moderately, but the longer one can acquire a huge range of values 
depending on the actual arrangement of the system as demonstrated in Fig. [5| 

The full picture of the excitation dissipation depends essentially on the interplay between the relaxation to the 
ground state and the thermalization driven by the system-bath interaction. The total process, which is a two- 
exponential decay, is shown in Fig. |6] The exponents can be attributed to the thermalization and relaxation to the 
ground state (cf. Eq. (26l). In the long-time limit, both the lower and the higher state populations decay with the 
same rate (the latter not shown here) because the equilibration is much faster than the relaxation. Therefore, the 
systems relaxes while being in a dynamic equilibrium. The reorganization energies partly determine the rates of 
the initial step of evolution and the amplitudes of the long-time decay. The long-time rates are identical for all four 
combinations of A's (as seen from the parallel slopes in the figure). These dependencies are qualitatively similar 
to those of the evolutions without relaxation to the ground state. Upon comparison of the numerical results with 
analytical expressions, it can be seen that the faster relaxation rate largely determines the rate of the whole process, 
and therefore the lifetime of the short-lived excitonic state is the determining factor of the excitation relaxation 
to the ground state. Moreover, as can be seen in Fig. [t] and as follows from Eqs. (27 1, in the presence of the 
system-bath interaction the decay rates are different depending on which of the two states is the lower one (cf. the 
case of purely excitonic mixing of lifetimes). The decay is faster if the short-lived state is below the long-lived one, 
however, for the parameters used, the difference is just of a few picoseconds. 

The parameters of the heterodimer chosen for calculations are similar to those of the Chl-Car dimer, which is 
assumed to be responsible for the additional quenching of the excitation under the NPQ conditions [B]. Indeed, the 
excitation lifetime of the Si state of carotenoids is short (less than 10 ps) and the optical transition from the ground 
state to the 5*1 state is forbidden, hence the excitation and reorganization energies in this case are not directly 
observable. However, the reorganization energy is a decisive parameter of the system-bath coupling, therefore, we 
have to consider the possible influence of the reorganization energy of the Car molecule on the excitation relaxation 
in such a heterodimer. It might be expected, that the excitation energy of the Car molecule in a specific Chl-Car 
dimer should be slightly above the excitation energy of the Chi molecule under unquenched conditions, and tends 
to be lower under the NPQ ones [8]. This statement can be examined using our model considerations, which 
demonstrate that it is only partially valid. The short lifetime component is the dominant control parameter of 
the excitation lifetime in such a dimer independently of the relative positioning of the molecular excitations of the 
monomers. As shown in Fig. [7] the excitation decays with the characteristic lifetime of 8ps if the optically allowed 
excited state of the Chi molecule is higher than the optically forbidden excited state of Car. It decays with the 
lifetime of 13 ps in the opposite case. Fast (of the order of 100 fs) relaxation of the excitation from the initially 
excited state of the Chi molecule to the excited state of Car taking place in the former case determines the main 
difference between both situations. Thus, the difference in the quenching efficiency between these two cases are 
related to the difference in the rates of excitation transfer back to the chlorophylls of the surrounding antenna. 
Moreover, the dominant effect of the quenching ability of such heterodimer is determined not by the excitonic 
effects (see Fig. [s]), but rather by the exciton-bath interaction. As can be concluded from Eqs. (27 1, the molecule 
characterized by the shortest relaxation time forms the valve of the excitation relaxation and defines the excitation 
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lifetime within the heterodimer. 

The same model of heterodimer is also applicable to exciton-CT state coupling conditions [211 1^ Since 
the lifetimes of both states are long and comparable, the dynamics within such a dimer corresponds to the limit of 
/C = 0. The quenching ability of this dimer is determined solely by the thermalization within the excitonic states. 
However, Fig. ^jp (or d) demonstrates that, under certain configuration of reorganization energies, the CT state can 
be the higher one (red curve) even though it should be the lower one to ensure the quenching ability in the Forster 
limit. 



5. Conclusions 

We have considered dynamics within a heterogeneous system of two interacting molecules. The main sources 
of heterogeneity in this system are the difference in excitation energies, asymmetry in reorganization energies 
(spectral line-widths) and the different excited state lifetimes of the chromophores. We point out that difference 
of reorganization energies introduces an ambiguity in the concept of energy gap between the excited states. This 
situation is particularly important when considering differences between coherent and incoherent transfer, since upon 
excitonic mixing under certain combinations of reorganization energies, the two states can become energetically 
swapped with respect to the uncoupled ones. Furthermore, employing the HQME technique revealed that the 
conventional excitonic basis becomes renormalized under the action of the bath. 

The analytical study of excitonic mixing of the excited state lifetimes and the estimates following from it revealed 
that the lifetime of the long-living state can be indeed greatly reduced. However, we found that the final picture 
of energy dissipation from the system crucially depends on the system-bath interaction rather than just on the 
excitonic coupling. Namely, the short-living state, irrespectively of its energetic position, forms the valve of the 
process and largely determines the rate of dissipation by setting its lower boundary. We conclude that in the case 
of the Chl-Car dimer the rate of excitation relaxation is of the similar order as the inverse lifetime of Car Si state 
even if the state is above the Chi excited state. This has direct consequences in assessing the model of Chl-Car 
dimer as a quenching center of the NPQ process. 
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Appendix A. Transformation of the Relaxation Superoperator 



The unitary transformation of a superoperator follows from the definitions of the analogous transformation of 
an operator and the action of a superoperator, accordingly: 



p^U ^pU Pab = '^{U ^)aiPijUjb] 

p = U^^p =4> Pab = y^^(M^^)ab,cdPcd- 



(A.l) 
(A.2) 



cd 



This way any superoperator O can by transformed as O = W ^OU. Upon comparison of (A.l I and (A.2) we can 
immediately deduce a rule for composing the unitary transformation superoperator: 
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(pi ^)ab,cd = {U ^pU)ab\p,j=5i,5ja- 

For instance, in the case of a heterodimer, we use Eq. (10 1 to obtain the following superoperator: 
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(A.4) 



here we use the shorthand notations: c — cos 6, s = sinO. Therefore, bearing in mind that ICab,ab — K,ba,ba — , 
it is straightforward to show that the result of a product U^^K.U for the population relaxation elements yields Eqs. 
pi. 
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